Disorder strength and field-driven ground state domain formation in artificial spin ice: 

experiment, simulation and theory 
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Quenched disorder affects how non-equilibrium systems respond to driving. In the context of 
artificial spin ice, an athermal system comprised of geometrically frustrated classical Ising spins 
with a two- fold degenerate ground state, we give experimental and numerical evidence of how such 
disorder washes out edge effects, and provide an estimate of disorder strength in the experimental 
system. We prove analytically that a sequence of applied fields with fixed amplitude is unable to 
drive the system to its ground state from a saturated state. These results should be relevant for 
other systems where disorder does not change the nature of the ground state. 

PACS numbers: 75.50.Lk, 75.10.Hk, 75.75.-c 



Artificial spin ice [1-^] consists of nanofabricated ar- 
rays of elongated magnetic dots that are small enough to 
be single-domain, but large enough to be athermal. Each 
dot can be represented by a macro Ising spin, whose dy- 
namics is governed by frustrated magnetostatic interac- 
tions and an external magnetic field. Like other artificial 
systems [^-^ij, artificial spin ice offers a setting in which 
to explore athermal, nonequilibrium dynamics, with the 
advantage that its microscopic configurations are easier 
to measure than those of, e.g, granular materials. 

Understanding these dynamics is difficult as the study 
of nonequilibrium systems lacks a key ingredient: the 
probability distribution of microstates, i.e. the Boltz- 
mann factor, which allows the determination of all rele- 
vant physical quantities Furthermore, athermal sys- 
tems are intrinsically strongly out of equilibrium and the 
relevant distribution function is not obtainable from per- 
turbations of the equilibrium distribution. Artificial spin 
ice is an experimental example of such an athermal sys- 
tem, with the additional ingredient of geometrical frus- 
tration. Furthermore, it has been shown jsl-flU that al- 
though temperature is not relevant, randomness does en- 
ter via quenched disorder, due to small unavoidable vari- 
ations during fabrication. 

Previous studies have quantified disorder strength in 
artificial spin ices ll-[l3[ flBj , but relatively little has 
been said about how it affects dynamics, especially in 
square ices. In particular, ideal square ices have a well- 
defined, two-fold degenerate ground state [s", Tel (GS; 
see Fig. Hl^a)), but this has proven unattainable in ex- 
perimental studies of field-driven demagnetization, which 
have yielded states with only short-range GS correla- 



tions [ill, [HI . Simulation studies of a nanopatterned 
superconductor "spin" ice suggests that disorder is par- 
tially responsible [9], but its full role has not yet been 
elucidated. 

In this Letter, we examine the extent to which disorder 
can disrupt ordering processes in artificial spin ice. Our 
argument is a specific example of a broader problem of 
how disorder affects access to a set of degenerate states, 
and is also applicable to, e.g, antiferromagnetically cou- 
pled magnetic dots with a distribution of switching bar- 
riers [19-21]. Our results complement previous studies of 
how disorder affects phase space 0, ^[ . 

We present experimental studies in which rotating field 
protocols are used to drive dynamics in a square ice, in 
the manner simulated in Ref. |2J. In those simulations, it 
was shown that - in ideal systems - a constant-amplitude 
rotating field can generate states with large domains of 
GS ordering via orderly invasion processes. In our exper- 
iments, GS domains are smaller than predicted, which 
new simulations show is attributable to quenched disor- 
der. We prove analytically that disorder blocks GS access 
in real systems for any field protocol with constant am- 
plitude, by forcing multiple GS domains to nucleate. 

An important related problem is that of how the basic 
"step" taken by the system as it moves through phase 
space affects the accessibility of states: for example, in 
Monte Carlo simulations of vertex models [25], single spin 
flip dynamics can "freeze" in regions of phase space far 
from the ground state, an issue not faced by loop-based 
dynamics. Indeed, this problem is general, appearing also 
in contexts such as domain wall creep dynamics [26]. In 
this present Letter, we address the question of pathways 
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FIG. 1. (a) A pure tiling of Type 1 vertices gives a ground 
state (GS) configuration. The other GS is obtained by flipping 
all spins, (b) One of four diagonally polarized states (DPS) 
and the field direction used to obtain it. (c) The 16 vertex 
configurations, grouped in order of increasing energy. 



to the GS in the context of single spin flip dynamics, an 
approach motivated by the stochastic nature of spin flips 
in artificial spin ice, but leave open the possibility that 
dynamics governed by global moves may be different. 

As in Ref. 24, we study two array edge geometries. In 
"open edge" ("closed edge") arrays, edge spins have an 
odd (even) number of nearest neighbors. These geome- 
tries are shown in the insets to Fig. [2fa,e). In perfect 
systems, differences in coupling (of the order of nearest- 
neighbor interactions) at array edges cause edge geome- 
try to "select" dynamics with distinct field dependence 
for the two array types. We will see that disorder disrupts 
this. 

Five nominally identical arrays of each edge type were 
patterned on a single Si chip with electron beam lithog- 
raphy, as per Ref. |28|. Islands were nominally 85 nm by 
280 nm on a lattice of 400 nm constant, with a thin film 
structure of Cr(2nm)/Permalloy(30nm)/Al(2nm), form- 
ing moments of ~ lO^/i^, with nearest neighbor coupling 
of ^ 10 Oe. As is common in demagnetization stud- 
a large in-plane field H ■ 
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29|, a large in-plane field H = 2 kOe was applied 
along a diagonal symmetry axis to prepare a diagonally 
polarized state (DPS; see Fig. [TJb)). The field was then 
reduced to a hold value i^h and the sample was rotated 
in-plane. We studied hold fields at 22 Oe increments 
between 411 Oe and 606 Oe. After hundreds of rota- 
tions - more than enough to reach a predicted steady 
state ^4] - the field was ramped to zero at a rate of ~ 
10,000 Oe/s (compared to a rotation period of ~ 30 ms). 
We have confirmed by simulation that this ramp-down 
does not cause the demagnetization effects seen for slow- 
ramp protocols [8, 30], because the field range over which 
non-trivial dynamics can occur, i^max — ^min ^ 150 Oe, 
is crossed within a single rotation. For large i^h, the 
field angle at which ramp-down starts influences the fi- 
nal configuration, as described below, but its exact value 



is unimportant in general. For each iJh, remanent states 
were imaged by MFM, which shows a pole at each is- 
land end, indicated by the red/blue contrast in Fig. [21 
confirming that the islands are single-domain. 

Array configurations can be conveniently represented 
in terms of vertices. Figure [TJc) shows the sixteen pos- 
sible vertex configurations grouped into types based on 
energy jl| . The GS is a chess-board tiling of Type 1 ver- 
tices, and the DPS that we use as an initial configuration 
is a uniform tiling of one of the Type 2 vertices. Dynam- 
ics on DPS or GS backgrounds can be described in terms 
of the motion of Type 3 vertices [Til, [H, [111, ill . The 
population of Type 1 vertices serves as an indicator of 
the level of GS ordering 0, 

Figures [2] (a-d) and (e-h) show example MFM images 
from the i^h series, for open and closed edge arrays re- 
spectively. Key configurations are mapped schematically 
(insets) in terms of dipole moments and vertex type. 
Note that images for different Hi^ are not all taken from 
the same array. Average fractional vertex populations 
for each i^h are tracked in Fig. [3fa). We make two gen- 
eral observations. First, the population statistics and net 
magnetizations (not shown) are almost identical for open 
and closed arrays for all i^h, indicating that edge effects 
are suppressed, which we show below is caused by disor- 
der. Second, the maximum Type 1 population, found at 
i^h = 520 Oe, is significantly reduced compared to the 
predictions of the ideal model [24] . 

Examining the configurations attained, we see that for 
i^h < ^min ^410 Oe, the field does not affect the initial 
DPS. For 410 Oe < i^h < 455 Oe the configurations 
at remanence consist of chains of reversed moments on a 
background of the initial DPS, as seen in Fig.[2fa,e). It is 
clear from the MFM that most chains are nucleated in the 
bulk, presumably at sites with low switching barrier 
The chains are similar to those reported previously in dc 
field experiments 1113, 15, 32|. In those experiments, 
they occur via bulk nucleation, cascading and pinning 
under the influence of interactions and disorder. 

As i^h is increased to 476 Oe, smah GS domains 
form, and mixed GS/DPS phases are found, as seen in 
Fig.[2](b,f). Near 520 Oe, increasing numbers of moments 
are reversed from the initial DPS, the magnetization ap- 
proaches 0, and all four Type 2 vertices reach similar pop- 
ulations - the memory of the initial DPS is lost. The GS 
domain size increases, and the Type 1 populations reach 
a maximum of 50%, as seen in Fig. [2jc,g) and Fig. [3fa). 
An example GS domain is outlined in Fig.[2fg). Domain 
wall structures separating GS domains, similar to those 
caused by thermal ordering [28], are also observed. 

Increasing i^h further rapidly suppresses GS order, as 
Zeeman energy dominates and DPS ordering that couples 
to the field is preferred. Only moments with large switch- 
ing barriers can pin; the rest align with the field. This 
is evident in the increasingly polarized states observed in 
Fig. [2l(d,h). The magnetization direction is determined 
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FIG. 2. MFM images of final states of open (a-d) and closed (e-h) arrays after rotation at selected hold fields. Insets indicate 
the corresponding island moments for selected areas of the array. A GS domain is outlined in (g). Images are false-colored 
using the software package WSxM [27l l. 
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FIG. 3. Vertex populations vs hold field for (a) experiment 
and (b) theory. Symbols represent vertex types as shown in 
the legend, with open (closed) symbols for open (closed) edge 
arrays. Each data point is the average over several runs; error 
bars represent the standard error. 



by the field direction at which ramp-down occurred. We 
estimate i^max to be 560 Oe; above this field, no sig- 
nificant GS domains are found, and simulations indicate 
that any GS ordering is picked up during ramp-down. 

We now use numerical simulations to establish that the 
above observations can be explained by quenched disor- 
der, and to estimate its strength relative to other energies 
in the system. In our simulations, the Ising spin i flips if 
the total field acting on it, comprising the external field 
and dipolar interactions with all other islands, exceeds 



(1) 



where rhi is a dimensionless unit vector along the spin di- 
rection and h^c^ is the island's switching barrier, which, 
in a perfect system, is the same for all islands. This 
threshold-based model [H, [ll|, [li', ^33] has a cos 6 angu- 
lar dependence, and is appropriate for describing the 
Zeeman-energy-driven propagation of domain walls, such 
as occurs during reversal of dots with dimensions similar 
to ours [34], in which domain wall nucleation is assisted 
by the curling of magnetization at island ends 12911 . (For 
smaller dots, S toner- Wohlfarth switching fsH*, 36] would 
be more realistic.) Like other authors [11-13, 15], we im- 
plement disorder by taking the h^c^ from a Gaussian dis- 
tribution with standard deviation a; we show elsewhere 
that this type of disorder behaves similarly to disorder 
in interactions [37]. We work in reduced units where 
the nearest-neighbor dipole coupling is 1.5 relative to 
M^/47r/io is the island net moment) and the mean 
he value is 11.25. 

We simulate a protocol in which the field rotates with 
constant amplitude h and angular step (i6> = 0.01 radians 
for 10 cycles, long enough to obtain a steady state. In 
line with experiments, the field is then ramped down over 
half a cycle to = 8, a field strength too low to induce 
dynamics. At each field application, the system evolves 
by flipping single spins according to criterion ([1]) until no 
further flips are possible. 
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We find good agreement with experimental vertex pop- 
ulations - in terms of general trends, peak rii value, and 
lack of dependence on edge geometry - when disorder is 
in the "strong disorder regime" of Ref. 37. For example, 
Fig.[3l^b) shows results for a = 1.875, a distribution width 
equal to 125% the nearest neighbor coupling, and large 
enough to suppress edge effects. This value is in agree- 
ment with the value of cr = 60 Oe, relative to a mean 
switching field of 320 Oe, given by Pollard et al [15], who 
studied arrays similar to ours. 

Disorder in our simulations is an effective switch- 
ing dispersion that incorporates effects from disorder in 
switching characteristics and interactions. As a point of 
comparison, if there was no disorder in interactions and 
island critical fields were directly proportional to their 
volume, a = 1.875 would correspond to a standard de- 
viation in island linear dimensions of 5%. Alternatively, 
if disorder originated only from fiuctuations in nearest- 
neighbor interactions, the disorder would correspond to 
a standard deviation of 40% [37]. We have been able 
to measure the standard deviation in island dimension: 
the value ~ 1% indicates that both types of disorder are 
present. 

Having seen that a rotating field does not drive a dis- 
ordered system from the DPS to a GS, we now ask: to 
what extent is this inherent to the nonequilibrium driven 
dynamics of a frustrated system, and what is the role 
played by disorder? Here we prove that when disorder 
is present no protocol with constant field amplitude can 
force a single GS ordering to cover the array, if disorder 
is strong enough or the system is large enough. 

A key ingredient of our argument is the two-fold degen- 
eracy of the GS, which allows for separate GS domains 
to form. Thus, we require that unlike in, e.g, the ran- 
dom field Ising model (in which the GS becomes more 
accessible for stronger disorder [38| ) , disorder should not 
change the nature of the GS. This is true for switching 
field disorder. 

We consider two mechanisms by which formation of 
separate GS domains can occur. We calculate an up- 
per bound on the probability P(not blocked) that neither 
mechanism operates; P (blocked) = 1 — P(not blocked) is 
a lower bound on the probability the GS is blocked. We 
outline the argument here, and give details as Supple- 
mental Material [1^. 

The first mechanism depends on the initial state being 
a DPS. To drive a DPS to a GS, half the spins must be 
fiipped, and the spins of the DPS can be divided into two 
groups based on their alignment with either GS. Suppose 
one spin from each of the two groups is pinned and re- 
mains always in its initial state, e.g, spin C and D in 
Fig. [TJb) . A single GS cannot contain both C and D 
in their initial states, but a configuration with two GS 
domains can. 

The GS is not blocked by pinning only if at least one 
of the two groups contains no pinned spins. Then, the 
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FIG. 4. The minimum GS blocking probability vs (a) 
disorder strength for arrays containing 400-1000 spins; and 
(b) array size for disorder strengths from 1.0 to 2.0. Inset: 
P (blocked) vs field strength, for a 20 x 20 array with a Gaus- 
sian distribution of he with standard deviation 1.875. Numer- 
ical values are in the same reduced units used in the rest of 
this Letter. 



first terms in the expression for P(not blocked) depend 
on the probability, Ppin, of a spin being pinned - that 
is, the probability the spin's switching barrier is so high 
it cannot flip, even in a maximally unfavourable local 
environment. This depends on field strength h. Two 
limits are Ppin ^ 1 as ^ and Ppin ^ as ^ oo. 

The second mechanism of GS blocking does not rely 
on an initial DPS. If two spins that are antiparallel in 
the GS - e.g, spin A and B in Fig. [TJa) - are "loose" 
and always align with the external field, then the GS is 
blocked. This gives a second set of terms in the expres- 
sion for P(not blocked). Pioose, the probability a spin 
aligns with an external field even when its neighbors are 
GS ordered, has the limits Pioose ^ 1 as ^ oo and 
^loose and h ^ 0. 

The probability that the GS is blocked is 

P(blocked) = 1 - [2(1 - Ppi„)"/2 _ (1 _ 

X [4(1 - Ploose)"/'(l - (1 - Pioose)"/')' 

+4(1 - Pioose)'"/'(l - (1 - Pioose)"/') + (1 " Pioose)"] • 

(2) 

The inset to Fig.SJ^a) shows P(blocked) vs field strength, 
for the system studied in simulations. P(blocked) > 65% 
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always. Figure |4] shows that the minimum of P(blocked) 
grows rapidly with disorder strength and array size. In 
the limit of an infinite system, finite probabilities of 
pinned and loose spins lead to finite populations of spins 
in both GS alignments, and the GS is necessarily blocked. 

Because we have been conservative in our estimates 



of 

Ppin and Piq 



these results are a lower bound on 



P(blocked). While large P(blocked) indicates the GS is 
inaccessible, small P(blocked) does not mean it can be 
reached: for example, ideal systems can jam [24]. Our 
results apply to any field protocol with fixed field ampli- 
tude, such as field protocols where the sense of rotation 
alternates. An open problem is whether protocols with 
varying field amplitude face similar blocking. Finally, 
we have shown that although rotating fields do not at- 
tain the GS, they do achieve a high level of GS ordering, 
pointing to questions about interplay between disorder 
and optimization [38i, iol, 41|. 
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